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Abstract. In this review presented at the Symposium A stellar journey in Uppsala, 
June 2008, I give my account of the historical development of the MARCS code from 
the first version published in 1975 and its premises to the 2008 grid. It is shown 
that the primary driver for the development team is the science that can be done 
with the models, and that they constantly strive to include the best possible physical 
data. A few preliminary comparisons of M star model spectra to spectrophotometric 
observations are presented. Particular results related to opacity effects are discussed. 
The size of errors in the spectral energy distribution (SED) and model thermal 
stratification are estimated for different densities of the wavelength sampling. The 
number of points used in the MARCS 2008 grid (108000) is large enough to ensure 
errors of only a few K in all models of the grid, except the optically very thin 
layers of metal-poor stars. Errors in SEDs may reach about 10% locally in the UV. 
The published sampled SEDs are thus appropriate to compute synthetic broad-band 
photometry, but higher resolution spectra will be computed in the near future and 
published as well on the MARCS site (marcs.astro.uu.se). Test model calculations 
with TiO line opacity accounted for in scattering show an important cooling of the 
upper atmospheric layers of red giants. Rough estimates of radiative and collisional 
time scales for electronic transitions of TiO indicate that scattering may well be the 
dominant mechanism in these lines. However models constructed with this hypothesis 
are incompatible with optical observations of TiO (Arcturus) or IR observations of OH 
(Betelgeuse), although they may succeed in explaining H 2 line observations. More 
work is needed in that direction. 



PACS numbers: 97.10.Ex , 97.20. Jg , 97.20.Li , 97.20.Pm , 97.20.Tr, 95.30.Ky) 
1. MARCS : little history 

The MARCS model atmosphere code has been in use since the mid-70's. It has its roots 
in Bengt Gustafsson's early work on a Feautrier-type method for model atmosphere 
calculation including convection (Gustafsson 1971). A number of features, including 
detailed continuous opacities were already included in that pre-MARCS version. And 
in a Gustafsson way it was directly applied to an astrophysical problem : determining 
the metallicity of F type stars (Gustafsson & Nissen 1972). A few years later the 
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collaboration with Eriksson, Nordlund, and Bell gave birth to MARCS: a code for 
Model Atmospheres in Radiative and Convective Scheme (Gustafsson et al. 1975, Bell 
et al. 1976). It allowed the computation of hydrostatic, plane-parallel (PP), line- 
blanketed atmospheres, with convection included following Henyey et al.'s (1965) recipe 
for MLT. Line opacity was included in the form of Opacity Distribution Functions 
(ODF), and FGK type stars, also metal-poor were modelled. The 1975 paper is a 
highly recommended reading in the field, as it discusses in great details the underlying 
assumptions, the algorithm, its implementation, and the resulting model structures 
with the impact of varying stellar parameter^} The grid was also used right away 
in astrophysical applications: abundance determinations, colour calibrations, ... This 
version of MARCS models have been used in almost 2 decades by the MARCS creators 
and collaborators, and probably longer by others. 

Simultaneously, there were incentive to (i) update the code and its input data, 
and (ii) extend its applicability to other spectral types. The latter became effective 
with Eriksson et al.'s (1984) efforts, who less than a decade later could deliver cool 
carbon star models with polyatomics included (HCN and C2H2), and apply them to the 
determination of abundances of carbon stars (Lambert et al. 1986). This development 
of MARCS was made possible because the expertise necessary to compute these new 
opacities was included into the team. Eriksson et al. (1984) could demonstrate the strong 
effect of polyatomics opacities on the atmospheric structures and emergent spectra, and 
already pointed out the necessity to refine these opacities [§] Two years later carbon 
star models were computed also with the Opacity Sampling method by Ekberg et al. 
(1986). They demonstrated that the ODF hypothesis does not work well when opacities 
from different sources and not correlated in wavelength (e.g. diatomics and polyatomics) 
dominate at various atmospheric depths. The OS method is then more reliable provided 
a large enough number of sampling points are used (see Section 4.1). At that time 
computer limitations hampered the wide use of OS, but the trend was launched: next 
MARCS would use OS. 

The extension towards cool oxygen-rich stars is where it all started for me. I had 
started my PhD in France on the empirical modelling of Mira stars based on speckle 
interferometric observations. As I was getting about nowhere, my supervisor had the 
good idea to send me to Uppsala, to spend about a year with Bengt Gustafsson and 
learn as much as possible about cool star atmospheres. The challenge became quickly 
to produce cool star atmospheres for oxygen-rich stars, as only a carbon star setup was 
available^ Opacities for TiO and H2O were missing, the chemical equilibrium did not 

% incidentally, I was 12 then, and as far as I can remember this is about the time I got interested in 
astronomy, but the reason must be in a, at least initially, more contemplative approach to the stars... 
§ Thanks to the hard work of our colleagues in quantum chemistry this problem is alleviated for HCN 
(Harris et al. 2002, Harris et al. 2006), but is unfortunately still pending today for C2H2. I again 
strongly encourage molecular physicists and quantum chemists to assemble a line list for that species. 
There is a great reward : analysis of carbon AGB stars in the local group and the history of carbon, 
understanding of cool carbon dwarfs, ... 

|| This turned out to demand much more time, and I stayed in Uppsala a total of 3.5 years: first for 16 
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include important species like TiO, and sphericity was included in a trial version but 
did not work properly. Thanks to Mats Larsson, Bosse Lindgren, and Lars Pettersson 
at Stockholm university I was introduced in the arcanes of quantum chemistry, and 
could produce a first line list for TiO. Together with John Brett we computed VO and 
CaH line lists. We could not compute H 2 in the same manner, as our limited quantum 
chemistry skills stopped at 2-atoms species, but used a trick to produce a pseudo line list 
from NASA observations of rocket exhausts (mean opacities, and line density, together 
with a distribution of strengths). We knew this was in need of improvement, as well 
as the other molecular opacity, but it would allow, just as for carbon stars, a first 
assessment of the impact of these species on the atmospheric structure. The inclusion 
of spherical symmetry was made with the help of Ake Nordlund using his algorithm 
(Nordlund 1984). Finally, The MARCS code was made cool oxygen-rich capable by the 
modification of the chemical equilibrium routine, and OS was substituted to ODF with 
11000 points, however keeping ODF as an option. The extensive Kurucz atomic line 
data was included by Bengt Edvardsson (see Edvardsson et al. 1993). In Plez et al. 
(1992), the improved code is described, as well as the impact of the molecular opacities, 
and of sphericity on models. Low resolution spectra are compared to observations of 
M giants and the models (Plez 1992) are compared to other grids. Again these models 
were used with success in many applications in the following years. 

Using the same code skeleton, Jorgensen et al. (1992) published a grid of carbon 
stars, with 5400 OS points, sphericity, and C3 opacity, although the latter would prove 
to be overestimated by a large factor later on. They also discuss in details the effect of 
sphericity. 

With the same code as Plez (1992), Edvardsson et al. (1993) computed a grid of 
metal-poor solar-type stars with a combination of 4100 OS points in the UV (A < 4500A) 
and 1400 ODF points in the red, and applied it to the chemical analysis of Galactic disk 
stars, in a milestone paper. 

The evolution continued towards more exotic spectral types with the first grid of 
line-blanketed H-deficient models by Asplund, Gustafsson, Kiselman & Eriksson (1997&) 
(see also Asplund, Gustafsson, Kiselman & Eriksson 1997a). Bound- free opacities were 
updated using data from the opacity project, and free-free opacities were added for 
carbon and Helium ions. These were of course used to analyse RCrB and related stars 
(e.g. Asplund, Gustafsson, Lambert & Kameswara Rao 1997). 

There are a number of other papers summarising other add-ons and updates : 
Brett (1995) published a grid of M dwarfs, Bessell et al. (1998 b) (see also Bessell 
et al. 1998 a), in their work on the calibration of Johnson-Cousins photometry used an 
updated MARCS grid of cool star models (both giants and dwarfs, including updated 
molecular opacities, a much more extensive chemical equilibrium, and more OS points. 
Plez et al. (2003) added ZrO opacities and computed the first S-type star atmosphere 

more months as a french civil servant part-time teaching French conversation to students and personnel 
of Uppsala university; then on a Swedish PhD fellowship. I met my wife during that time, and since 
then Sweden has become a second home-country to me. 
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grid with detailed blanketing (100000 OS points), and made a first attempt at classifying 
S stars using the synthetic spectra and colours, Gustafsson et al. (2003) published on 
the web an extensive updated grid of FGK type-stars down to very low metallicities. 
The two latter grids used a completely revised chemical equilibrium with 92 elements 
and their ions as well as over 500 molecules, with updated partition functions and 
dissociation energies. 

2. MARCS 2008 

The papers above, starting around 1998, were all paving the way, and announcing the 
much awaited new MARCSj^j The idea of a large grid of updated MARCS models 
covering the cool part of the HR diagram dates back to a suggestion Bengt made to 
me in the park on Rackarberget at Dan Kiselman's PhD party in 1993. Updates were 
constantly added, but culminated in the 2-3 years before 1998, while I was on the payroll 
of Uppsala Observatory again. Citations of Gustafsson et al. in preparation date about 
back to 1998. At some point we had a single version of the code, and a single set of 
physical input data, but I got my job in Montpellier and the entropy started increasing 
again. This is without mentioning the Copenhagen branch that diverged already in 
1991. We had a number of small workshops or get together along the years to try to 
(i) get a single version of the code, (ii) update input data, (iii) find tricks to avoid 
crashes of extreme models like high luminosity, low gravity, cool supergiants, (iv) draft 
a series of papers and decide what to discuss in what paper. Then of course in between 
these meetings, other exciting projects, teaching and administration took over. During 
the past 9 months, my sabbatical in Uppsala allowed us once more, and this time for 
good, to get a single more debugged, more updated (esp. continuous and line opacities) 
code from the Montpellier and Uppsala ones. The first paper in a presumably long 
series appeared recently (Gustafsson et al. 2008) and the grid is (partly) on the web 
at marcs.astro.uu.se. MARCS 2008 is characterised by new opacities for H2O, atomic 
collisional line broadening included using the description of Anstee & O'Mara (1995), 
and hydrogen lines modelled using a code by Barklem, described in Barklem & Piskunov 
(2003). About 108000 OS points are used. Full details are provided in the paper, 
that also relates some of the historical background, and discusses in depth the physical 
assumptions, numerical methods, and physical data used. About 10000 models were 
computed 

% it was dubbed MARCS35 (guess why!) in 1998, and was supposed to be the final version for the new 
generation MARCS grid while we were sketching the first draft of a series of papers, the first of which 
appeared only 10 years later... 

+ Thanks to this care for details and documentation that characterises Bengt's work, I could compare 
computing (CPU) times for the pre-MARCS code in 1972 : 25mn for a PP model with 148 wavelengths 
accounting for the opacity of 25 Balmer lines, and MARCS 2008: lOmn for a spherical model with 
108000 points accounting for over 10 s lines. The computing time is comparable, but the physics is 
much more detailed! 
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3. Comparisons to observations 

One decisive test to be passed by the models is a comparison to observed high-resolution 
spectra and spectral energy distributions (SEDs). High- resolution spectra allow to check 
in detail if individual lines are well accounted for in the model, whereas low-resolution 
SED comparisons permit to assess more easily problems in continuous opacity, thermal 
gradients, and missing opacities. These comparisons must be carried out using stars 
with well determined stellar parameters spanning the grid of models. In this volume, 
Bengt Edvardsson presents the first comparisons to the Sun and solar-type stars SEDs. 
I have made a small number of comparisons for M-type dwarfs and giants. The results 
I present here are very preliminary. The model SEDs are the fluxes sampled at the 
108000 OS points and as such are not like an observed spectrum (see Figure Q. When 
averaging on a wavelength interval including many OS points, the mean average level 
of the spectrum is recovered as, statistically, the OS points hit the continuum and lines 
of all strengths. So, the sampled fluxes are a good representation of the SED only when 
binned to a much lower resolution. Our OS scheme has a constant wavelength resolution 
of A/AA = 20000 between 900A and 20/im. 

Comparisons were made with spectra from the MILES (Sanchez-Blazquez et al. 
2006) and STELIB (Le Borgne et al. 2003) libraries, with both observations and 
synthetic SEDs binned to 50A intervals. For the giants (Ml to M6III, i.e. 3800K 
to 3250K), the fit is generally good, except for excess absorption in the models around 
5100A most probably due to MgH, and ups and downs shortward of 4200A reminiscent 



of what Edvardsson (this volume) finds for the Sun (see Section 4.1 for a possible 
explanation of part of this). At M6III, there is a hint of a missing opacity around 
7500A (LaO?). The dwarfs are more difficult. Using a T e ff=3950K allows a fit of the 
TiO bands of an M0.5V template, but the computed spectrum is then too blue. A lower 
temperature of 3750K gives a better shape for the continuum, at the expense of too 
deep TiO bands. In all cases there is excess blue flux in the model. At M6V, TiO bands 
are fitted but there is also excess flux in the blue. Work with Mike Bessell in the few 
days after this conference, using his own higher resolution spectrophotometry, allowed 
us to tie the blue flux problem to CaOH absorption, and maybe A1H in metal-deficient 
stars, while providing better global fits. Detailed spectral comparisons need to be done, 
especially for dwarfs, missing opacities should be assessed and if possible problem cured, 
and a temperature scale based on the models should be derived, as has been done for 
red supergiants (Levesque et al. 2005, Levesque et al. 2006). A paper in the MARCS 
2008 series will be devoted to M stars and will detail these comparisons. 



4. Opacity effects 

A thorough discussion of the general effects of blanketing, of the impact of different 
opacity sources, and of the effect of abundance and microturbulence changes on 
blanketing in various temperature regimes of the MARCS models is presented in 
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Figure 1. Illustration of the fact that a sampled SED, as coming out from the 
MARCS code, is not equivalent to a spectrum at the same resolution. The MARCS 
flux is computed exactly at a number of prescribed wavelengths, but nothing is known 
about what happens between the points. Nevertheless when binned at sufficiently low 
resolution, the sampled SED tends towards the spectrum. 



Gustafsson et al. (2008). Here I will only complete this discussion with two items: the 
impact of wavelength sampling, and some trial calculations on scattering in molecular 
lines. 

4-1. Wavelength sampling 

I already stressed the fact that a large number of OS points is necessary in order to 
statistically well represent the opacity, i.e. to sample equally well the continuum, and 
lines of all strengths. But what is large enough? There are two issues: (i) how many 
points are needed to give a good representation of blanketing effects, and therefore 
provide a converged temperature stratification, in the sense that more OS points will 
not lead to changes in the model stratification larger than a prescribed value (say IK), 
(ii) how many points are needed to give a good representation of the SED, with a given 
model stratification (Figure [l]). These two numbers may well be different. I computed 
series of models with subsets of our 108000 OS-points: 3 models with about 36000 
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points, 10 with 11000 points, and 30 with 3600 points. This was done for (i) a model 
of the Sun, (ii) a 3500K M giant, and (iii) a 5500K metal-poor dwarf ([Fe/H]=-2). For 
each set of models, the average spectrum and standard deviations were computed, and 
the same was done with the thermal structure (T — tr oss ). 

In the case of the Sun, temperature fluctuations are less than 10K above tr oss =- 
2.0, and less than 1 or 2K below, with 36000 OS points. They are of the order of 
30-50K above 7r oss =-3.5, with 11000 points (up to 100K for the most deviant model). 
Fluctuations stay below 15K in all cases for tr OS s > —2.0. The situation is much better 
for the cool giant, with temperature errors less than 20K everywhere (—6 < tr oss < +2, 
down to 3600 sampling points. Errors decrease to less than 4K everywhere with 36000 
points. On the contrary, the situation is worse in the optically thin regions for the 
metal-poor dwarf: temperature errors increase from 15K at Tr oss =-4 to 120K at r Ross =- 
6, even with 36000 points. The situation is better in the line forming region: less than 
10K fluctuations below r Ross =-2.0, at all samplings. So, errors in the T-structure due 
to sampling, are only a few K in all cases, for sampling densities in excess of 36000 
points. At smaller optical depths, errors may be large for metal-poor stars, and are 
small for cool giants (<10K). Errors in the temperature stratification occur only if the 
line opacity is not well sampled, i.e. if the line and the OS-point densities are both too 
low. 

Errors in the fluxes are below a few percents everywhere the flux matters when 
using 36000 points, with the notable exceptions of the UV flux of solar-type stars and 
the IR CO bands of red giants. The case of the Sun is illustrated in Figure [2} The errors 
around 3000A are of the order of 5 to 20%, and around 5% at 4000A. The trend at all 
wavelengths when increasing from 3600 points to 36000 indicates that the new MARCS 
grid with 108000 sampling points should provide sampled SEDs with systematic errors 
a/3 smaller. The cool giant model sampled SED shows errors of 5 to 10% in the IR CO 
bands at 1.6 and 2.5/iin, with 36000 OS-points. This is due to the intense CO lines that 
are sparsely distributed, and sampled at a resolution of about 6500. Not surprisingly, 
the flux of the metal-poor dwarf is very well modelled even with a low sampling (there 
are almost no lines), except where there are lines, i.e. in the UV. At the 36000 points 
sampling, errors reach more than 10% only below about 1600A, and are less than 2 to 
5% between 3000 and 4000A. 

In conclusion, the 108000 OS points used in the computation of the MARCS 2008 
models are sufficient for getting the errors on the thermal stratification down to a few 
K for all model parameters of the grid. An exception is the optically very thin layers of 
metal-poor stars, due to the fact that only a few lines dictate the thermal equilibrium 
in these layers. This is however a region where NLTE or 3D effects most probably play 
a greater role. The situation is worse for the sampled SEDs that may still be off by over 
10%, esp. in the blue-UV, as many lines act on the spectrum without greatly affecting 
the temperature structure, e.g., atomic lines in cool giants, see Fig. 5 in (Gustafsson 
et al. 2008), MgH or CH in the Sun. Of course better SEDs can be computed afterwards 
at very high resolution, and this will be done for the MARCS models. Sampled SEDs 
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Figure 2. Errors in sampled fluxes for Solar models computed with different numbers 
of wavelength points. The figure shows the standard deviation to the mean value of 
the flux in bins of AAA = 200 for each series of models with 3600, 11000, or 36000 OS 
points. The reference MARCS 2008 models are computed with 108000 OS points. 



can nevertheless be used, e.g., to compute broad band colours. 



4-2. Scattering versus absorption 

Molecular lines have a large effect on the thermal structure of cool stars, through 
heating or cooling of the outer layers, and backwarming of deeper layers. TiO numerous 
electronic transition lines in the optical are well known to lead to a large heating of 
the outer layers of cool stars. This is due to this strong opacity appearing on the blue 
side of the local source function and affecting thermal heating: q = J °° k x (J\ — B\)d\, 
(with q = in LTE, hot radiation from below, J\, and the large k\ in the blue impose 
an increase of B\ in the surface layers). Surface heating or cooling only happens if the 
opacity is in absorption, and a scattering term has no effect. It is well possible that 
some molecular lines do indeed form closer to scattering than pure absorption in the 
tenuous outer layers of red giants. Hinkle & Lambert (1975) already suggested, based 
on meagre available collisional excitation data, that electronic transitions of molecules 
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in cool star atmospheres may be dominated by radiative processes, and not by collisions 
with electrons or hydrogen. We may try to examine again the particular situation 
of TiO. Radiative transitions occur with t ra d = 40 to 100 ns in most TiO electronic 
transitions. This means that collisional rates of excitation between these levels must be 
at least 2.5 x 10 7 s _1 in order to compete and ensure LTE populations. For collisions 
with electrons if we use the, most probably inappropriate but only available we have, 
approximations of van Regemorter (1962) and Jefferies (1968) (see also Rutten 2003), an 
electron density of at least 3 x 10 14 cm~ 3 is needed. In typical MARCS cool supergiant 
models 2 x 10 6 < N e < 2 x 10 10 cm~ 3 , well below the limit. For hydrogen, Lambert 
(1993) proposed a modification of Drawin's (1969) formula, for hydrogen collisions with 
atoms that we may use, in lack of better estimates, to derive a critical hydrogen 
density of about 2 x 10 22 m~ 3 . In the same MARCS models as above, we find 
2 x 10 17 < Nh < 2 x 10 21 m~ 3 . Alternatively, using measured quenching rates in 
oxides, e.g. the measurements of Badie, J. M. et al. (2008) in YO, for Ar, He and O2, 
we find a critical density of perturbers of about 10 25 m~ 3 in order to reach collisional 
excitation rates of 2.5 x 10 7 s" 1 . The efficiency of H collisions would have to be at 
least 5 orders of magnitude larger than that of He for hydrogen to be of significance in 
the excitation of YO electronic levels, and similarly of TiO. It seems thus very likely 
that radiative processes indeed dominate over collisional processes in the population of 
electronic levels of TiO and in transitions between them. 

I therefore conducted a test calculation to study the impact on the thermal structure 
of red giant models if TiO lines were formed in scattering. The first set of models had 
T e ff = 4300 K, logg = 1.75 [Fe/H]=-0.5, and M=l M representing Arcturus. The 
model with TiO in scattering has a temperature about 250 K lower at r 500 o = 10~ 5 , the 
difference decreasing to about 130 K at r 50 oo = 10~ 4 and vanishing at r 50 oo = 10~ 3 , the 
pressure stratification remaining unchanged. This cooling is what Ryde et al. (2002) 
propose to explain the appearance of H 2 lines at 12/im in Arcturus spectrum. However, 
when using the cooled MARCS model to generate a synthetic spectrum in the optical, 
TiO band heads (e.g. 7(0, 0) at 7054A) b ecome visible, whereas they remain undetected 
in Hinkle et al.'s (2000) FTS spectrum of Arcturus. I did the same experiment with a 
model appropriate for Betelgeuse (T e g = 3600 K, logg = 0.0, [Fe/H]=0, and M=15 M . 
In this case the cooling amounts to about 250K at r 50 oo = 10 -5 , to about 130 K at 
r 5000 — 10~ 3 , and vanishes around r 50 oo — 2.5 x 10~ 2 . Ryde et al. (2006) detected 
the 12/xm water vapour lines in Betelgeuse as well but could not reproduce them with 
a standard atmosphere model at 3600 K, which is the temperature derived from the 
optical spectrum. They needed a much cooler (3250 K) atmosphere. Our cooled MARCS 
model does not, however, reproduce the water lines : they remain too faint. In addition 
neighbouring OH lines become too strong, but that could be alleviated by a decrease 
of O abundance. The 2.3/im CO bands (Wallace & Hinkle 1996) are marginally better 
reproduced with the cooled MARCS model, but not better than using a molsphere, i.e. a 
spherical shell of gas at 2000 K, on top of the photospheric model, as has been advocated 
by various authors to explain spectroscopic and high angular resolution observations of 
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red supergiants in the IR (e.g. Tsuji 2000, Perrin et al. 2007). 

Despite collisional time scale estimates being too long compared to radiative time 
scales to ensure LTE in electronic transitions of TiO, the simple assumption of its 
transitions occurring in scattering is not compatible with existing observations of red 
giants and supergiants. In particular it cannot solve the 12/im H 2 lines puzzle. A better 
assessment of collisional rates in TiO, as well as a full NLTE treatment of the electronic 
transitions, taking into account optical depth effects, as the lines may become optically 
very thick, would be of great value. Also the coupling with hydrodynamics should be 
studied in detail as the balance between expansion cooling and radiative heating may 
be shifted by large amounts in the upper atmospheric layers (see the contribution by 
Wolfgang Hayek in this volume) 

5. Final thoughts 

After over 30 years, MARCS is now as mature as it can be within its approximations 
(hydrostatic, LTE, MLT). A much updated grid covering the red part of the HR diagram 
is now being published on the web (marcs.astro.uu.se), and the first paper in a series 
that will e.g. detail the behaviour of all models, compare the synthetic spectra to 
observations for a selected sample of template stars, or discuss synthetic photometry, 
has just appeared (Gustafsson et al. 2008). Although we are more and more turning our 
development activities towards relaxing some of the simplifying hypotheses upon which 
MARCS rests, it is likely that updates of the MARCS grid will be made available in the 
future, specially when better or additional opacities can be included. The drive behind 
the development of MARCS has always been the science that can be done with it, and 
that is unlikely to change. It is a tool that despite eluding the dynamical side of reality, 
includes very detailed accounts of opacities, and is flexible and light enough that it 
can be used to test many ideas (effects of line scattering on atmospheric stratification, 
impact of NLTE over-ionisation on H _ opacity, etc). In that respect it will remain 
irreplaceable for many years to come. 
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Discussion 

Q: (T. Lynas-Gray) I note with interest the difficulty in matching energy distributions 
for M dwarfs. Did you use the water line list by Schwenke and Partridge, or the 
one by Barber and Tennyson? 

A: The one by Barber and Tennyson. 



